Avalanches in Strained Amorphous Solids: Does Inertia Destroy Critical Behavior? 
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Simulations are used to determine the effect of inertia on athermal shear of a two-dimensional bi- 
nary Lennard- Jones glass. In the quasistatic limit, shear occurs through a series of rapid avalanches. 
The distribution of avalanches is analyzed using finite-size scaling with thousands to millions of 
particles. Inertia takes the system to a new underdamped universality class rather than driving 
the system away from criticality. Scaling exponents are determined for the underdamped and over- 
damped limits and a critical damping that separates the two regimes. Systems are in the overdamped 
universality class even when most vibrational modes are underdamped. 

PACS numbers: 45.70.Ht, 61.43. Bn 



Many slowly driven physical systems exhibit long qui- 
escent periods punctuated by rapid avalanches [l[. Phe- 
nomena as diverse as earthquakes, Barkhausen noise in 
magnetic materials, dislocation cascades in single crystal 
microcompression and fluid interface depinning [2|-|8[ dis- 
play power law avalanche statistics in seismicity, acous- 
tic emission, slip, stress drop or interface advance. These 
power laws reflect a non-equilibrium critical transition at 
the onset of motion. 

This Letter addresses a fundamental question about 
the effect of inertia on such critical behavior. Power 
law scaling has normally been observed in overdamped 
systems. Studies of underdamped systems suggest that 
any inertia may drive the system away from the criti- 
cal point 0, [To| . In sandpiles, the onset of motion ap- 
pears to become a hysteretic first-order transition |ll| . In 
the Burridge-Knopoff model, inertia leads to a growing 
importance of non-critical, system spanning events 12 1. 



These conclusions about the effect of inertia seem at odds 
with the observation of power law scaling in earthquakes 
and laboratory compression tests 0, where seismic 
waves and acoustic emission indicate that the systems 
are underdamped. 

Here, quasistatic simulations of sheared glassy solids 
are performed over a full range of damping rates. The 
results reveal a rich phase diagram. Different universality 
classes describe the overdamped and underdamped lim- 
its, but both are described by critical finite-size scaling 
relations. The transition between the two limits occurs 
at a fixed damping rate that appears to have its own 
scaling behavior. Overdamped scaling extends to surpris- 
ingly small damping rates, where nearly all vibrational 
modes are underdamped. The power law describing un- 
derdamped avalanches is close to the Gutenberg-Richter 
law and an excess of large events is observed that is sim- 
ilar to observations of individual fault systems. 

Since we are interested in the general question of 
how inertia affects critical behavior, we consider a two- 
dimensional binary mixture of particles that has been 
widely studied as a model amorphous system I3-15|. 



The particles may represent atoms, grains, bubbles, col- 
loids or volume elements of a deforming fault zone. As 
particle size increases, temperature becomes less rele- 
vant. We focus on the athermal limit because it allows 
clear identification of small avalanches, and because other 
work indicates that temperature may also drive systems 
away from criticality [la, E3] • 

Particles interact via the Lennard- Jones (LJ) poten- 
tial, U(r) — 4«o[(<Zy /r) 6 — (a,y/r) 12 ] where r is the 
magnitude of the vector r between two particles and 
the species i,j are of two types, A and B. The particles 
have diameter claa = 5/3 cibb = a and clab — 4/5 a. 
The LJ energy and force are taken smoothly to zero at 
r c = l.ba.ij using a polynomial fit starting at 1.2a^ [l8| . 
Both particle types have mass m and the number ratio 
Na/N b = (1 + V5)/4. The depth of the inter-atomic 
potential uo sets the energy scale of particle interactions. 
The natural unit of time is r = y ma 2 /uq. 

Initial states are prepared as in Ref. [18], but, as there, 
the protocol has little effect on steady state shear. After 
annealing, the system contains iV(~ 10 3 — 10 6 ) particles 
in a square unit cell with edge L — 27a to 875a. The den- 
sity p — 1.38a~ 2 and the pressure is near zero. A pure 
shear strain is applied to the system by changing the peri- 
odic boundaries while holding area constant. The strain 
rate |e| < 10 _6 r _1 between avalanches was adjusted for 
each L to ensure simulations were in the quasistatic limit 
where results depend only on strain interval and not in- 
dependently on time. When an increase in kinetic en- 
ergy indicated the onset of plastic deformation, |e| was 
reduced to zero to allow the avalanche to evolve with- 
out external perturbation. Shearing resumed after the 
kinetic energy dropped below 1% of of the background 
value during shear. 

In order to model the athermal limit, the kinetic en- 
ergy released during avalanches must be removed by some 
damping mechanism. Unless noted, a viscous drag force 
was applied to each particle F<j rog = — Tmv where v 
is the non-affine velocity. As the dissipation rate T de- 
creases, the dynamics changes from overdamped to un- 
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derdamped (inertial) dynamics. We also show results for 
energy minimization dynamics that correspond to T — > 
oo. Vibrational modes with frequency uj > T are under- 
damped, and it is useful to compare lo to the root mean 
squared or Einstein frequency uie = y (w 2 ) = 17r . 

We focus on the steady state achieved after plastic re- 
arrangements have erased memory of the initial state 
(e > 7%). While the system is trapped in a local en- 
ergy minimum, work done by the applied strain leads to 
a nearly linear rise in potential energy and shear stress 
a s . When the minimum becomes unstable, there is a 
rapid avalanche of activity that leads to a sharp drop by 
E in energy and by A<r s in shear stress. In the over- 
damped limit, the system is trapped in the next local 
energy minimum. When damping is reduced, inertia can 
carry particles over subsequent energy barriers to reach 
lower energy states. One dramatic consequence is that 
the mean energy sampled by systems decreases by 30% as 
r decreases. Indeed, there is almost no overlap between 
the ranges of energy sampled for the three damping rates 
studied in detail below, Ft = 1, 0.1 and 0.001. 

To quantify the distributions of energy drops, we de- 
fine the event rate R{E, L) as the number of events of 
energy E per unit energy and unit strain in a system 
of size L d with dimension d — 2. A sum rule relates 
R(E, L) to the distribution of stress drops if Aer s is 
much smaller than the mean shear stress (a s ) and the 
shear modulus fi is nearly constant (l9j . The energy 
dissipated by avalanches comes from work on the sys- 
tem during stress increases. We introduce an elastic en- 
ergy S = Aa s L d (er s ) /4fi associated with a given stress 
change, where (a s ) /4/i 0.02 for all T and is indepen- 
dent of L. In steady state, the sum over stress drops 
equals the sum over stress increases and energy conserva- 
tion requires: / dEER(E 1 L) = J dSSR(S,L) = L d (a s ), 
where the last quantity is the total work per unit strain. 
Given the scaling exponents obtained below, both inte- 
grals are dominated by large events and thus R(E, L) 
and R(S, L) must exhibit the same scaling for large 
avalanches. Direct comparison of E and S for individ- 
ual events shows that they are not proportional for small 
events and we consider the scaling of both distributions 
below. 

Avalanche distributions are analyzed with finite-size 
scaling methods [20j j that assume the maximum size of 
events ~ L a is limited only by the system size L. The 
equations are developed for E but also apply to S. The 
distribution is assumed to obey the finite-size scaling 
ansatz: 



R(E,L) — L^g(E/L a ) 



(1) 



Equation Q] produces power law scaling for E « L a 
if g(x) oc x~ T for x « 1. One finds R(E, L) oc V E~ T 
with 7 = /3 + olt. One might expect that the probability 
of a small event in a given region would be independent 
of system size. This would imply a hyperscaling relation 
7 = j3 + ctr = d and r = 2, given the above relation 
/3 + 2a = d. As shown below, r < 2 and the number 
of small events rises much less rapidly than L d in our 
simulations (7 < d). Thus as L increases large events 
suppress small events either by changing the local con- 
figurations so small events are less likely to occur, or by 
increasing the probability that the same local configura- 
tion will produce a large avalanche. To our knowledge, 
the same behavior is not observed in other systems that 
display power law avalanche distributions. For example, 
the probability of small events is proportional to the size 
of the interface in models of fluid invasion or domain wall 
motion 




where g is an unknown scaling function. Using the sum 
rule in the previous paragraph, one finds a scaling rela- 
tion (3 + 2a = d = 2 as long as the integral of xg(x) is 
well defined. 



FIG. 1: Size of stress drops normalized by L a plotted against 
(a) Yt with a — 0.85 and (b) TL/c with a = 1.6. Statistical 
errors are smaller than the symbols for the indicated L/a. 

We first show how T influences the scaling of large 
events with L. The mean size (S) is not well-defined be- 
cause of the diverging number of small events for r > 1. 
Instead we show the ratio (S 2 ) / (S) L a , which gives sim- 
ilar event sizes and scaling as other moment ratios. As 
shown in Fig. [lja), a = 0.85 collapses results for all L 
onto a universal curve for Ft > 0.1. Results for Ft > 1 
and energy minimization are statistically indistinguish- 
able. As Tt decreases from 1 to 0.1, the mean avalanche 
size increases in the same way for all L. Since loet = 17, 
systems remain in the overdamped universality class even 
when almost all vibrational modes are underdamped. 
The key factor is not whether modes are overdamped 
but whether inertia can carry the system over the next 
energy barrier in the energy landscape. A small V can 
trap the system in the nearest minimum if the difference 
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in successive energy barriers is small and/or the path in 
phase space to the next barrier is complicated. 

As Ft decreases below 0.1 in Fig. QJa), results for 
different L separate. As shown in Fig. [IJb), results in 
this underdamped limit can be collapsed using a = 1.6 
and scaling T by the time L/c for sound propagation 
across the system at shear velocity c = 3.4a/r. Data 
for L > 27a follow a common scaling curve until Tt 
exceeds 0.1, and the critical behavior changes to over- 
damped. For TL/c « 1, all sound waves in the system 
are underdamped and all results tend to the same value of 
(S 2 ) I (S) L a . As TL/c increases past unity, the longest 
wavelength modes begin to be overdamped. The scaling 
exponent does not change, but (S 2 ) / (S) L a decreases. 
This is consistent with avalanches being cutoff by the 
wavelength of the largest underdamped mode rather than 
the system size. For the damping used here, this wave- 
length scales like 1/T and any finite damping takes the 
system away from the underdamped critical point. How- 
ever, long wavelength modes are always underdamped as 
L — s- in disordered solids because dissipation mecha- 
nisms must be Galilean-invariant and only damp relative 
velocities |2l| • Weakly damped simulations with two dif- 
ferent Galilean- invariant thermostats 3, 22 1 gave statis- 
tically identical results to those shown in Figs. 1-3 for 
TL/c < 0.1. 

We now examine the avalanche scaling in more detail 
in the overdamped (Tt — 1) and underdamped (Tt — 
0.001) limits and at the crossover between them (Tt = 
0.1). Figure [2] shows R(S, L)/ L" 1 vs. S in an overdamped 
system (lower curves) . Power law scaling is observed from 
S ~ 0.2uo to a cutoff that increases with L. Data in this 
scaling regime are collapsed with 7 = 1.3 ±0.1. As noted 
above, 7 < d implies that the number of events per unit 
area at a given S is strongly suppressed as L increases. 

Figure [3] shows finite-size scaling collapses of over- 
damped data for both E and S (middle curves). The 
scaling factor in S was chosen so that the two quanti- 
ties are comparable for large events, but the correlation 
between E and S breaks down below 0.2uo- Both quan- 
tities are well described by common scaling exponents 
Table (JTJ) from 0.2wo up to the largest event sizes. 

At lower energies, R(S, L) saturates while R(E, L) fol- 
lows a different power law that changes slightly with L. 
It is easy to confuse this power law with critical scal- 
ing if one only has results for R(E, L) at small L. The 
noncritical power law dominates R(E, L) for the smaller 
system sizes L < 109a used in previous studies. This ex- 
plains deviations in the reported values of a 13, 1|1] and 
why some papers concluded there was no critical behav- 
ior 23|, 24 1 . Dahmen has recently suggested that r has a 
mean-field value of 1.5 for all d in overdamped systems 
[25] ]. Our result of 1.2 is lower, but substantially higher 
than the values of t < 1 that would be inferred from the 
noncritical power law region in R(E,L) [nm. 

Figure also shows R(S, L)/L 1 for underdamped sys- 
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FIG. 2: Scaling of R(S, L) with S for overdamped r = 1 sys- 
tems with 7 = 1.3 (lower curves) and underdamped T = 0.001 
systems with 7 = 1.2 (upper curves). Underdamped curves 
are multiplied by 100 to prevent overlap, symbols indicate 
L/a and symbol size is larger than statistical errors. 
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FIG. 3: (a) Finite-size scaling collapse of R(S, L) and R(E, L) 
for overdamped systems (Tt = 1), underdamped systems 
(Tt — 0.001) and critically damped systems (Ft = 0.1) using 
the exponents in Table HI Statistical errors are smaller than 
the symbols and successive curves are shifted up by 2 or 3 
decades to prevent overlap. 



terns (upper curves). As expected, inertia leads to much 
larger avalanches. Results in the scaling regime for 
S > 0.2wo collapse with 7 = 1.2. The distributions all 
show a power law decay followed by a plateau that moves 
to larger S as L increases. While this is different from 
the sharp cutoff in the underdamped case, the form of 
the scaling function g(x) in Eq. [T] need not be simple. 
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r 


T 


Q 


P 


7 


1.0 
0.1 
0.001 


1.2 ±0.05 
1.0 ±0.05 
1.5 ±0.1 


0.85 ±0.05 
0.85 ±0.05 
1.6 ±0.1 


0.3 ±0.05 
0.4 ±0.05 
-1.2 ±0.1 


1.3 ±0.05 
1.2 ±0.05 
1.2 ±0.1 



TABLE I: Scaling exponents determined for overdamped 
(r = 1) and underdamped (F = 0.001) limits and at the 
crossover between them V — 0.1. Quoted values satisfy the 
scaling relations f3 ± 2a — d and 7 = /3 + olt and errorbars 
are estimated from the quality of finite-size scaling collapses 
for E and S using other V and moments. 



Figure [3] shows finite-size scaling collapses for S and E 
with the same scaling exponents (lower curves). In both 
cases, results for large events from different L collapse 
onto a universal curve. There is a plateau over a fixed 
range of a little under a decade followed by a very rapid 
decrease. 

These results clearly imply that inertia does not de- 
stroy critical behavior, but does lead to a different uni- 
versality class. Results for Galilean invariant thermostats 
with weak damping show the same scaling behavior. 
While our system lacks the complexity found in earth- 
quake faults, it is interesting to note that r is close to 
the value of ~ 1.6 for the Gutenberg-Richter law [2rjj |. 
In addition, the distribution of earthquakes for a given 
fault system typically has an excess of large events that 
is similar to the plateau seen in Fig. [3] 26]. 

The final example we consider is the intermediate case 
of Tt = 0.1 that seems to represent a crossover between 
overdamped and underdamped scaling in Fig. [TJ Fig- 
ure [3] shows a finite-size scaling collapse of R(E, L) and 
R(S, L) (top curves). The results follow a power law with 
t = 1 over 6 decades or more. Similar scaling was found 
for intermediate damping with Galilean invariant ther- 
mostats, for different interaction potentials, for simple 
shear, and in preliminary studies of 3D systems. This 
suggests that Tt = 0.1 represents a multicritical point 
separating regions that flow to underdamped and over- 
damped fixed points. 

In conclusion, introducing inertia does not destroy crit- 
ical scaling of avalanches in quasistatic shear of dis- 
ordered solids. Systems continue to be in the over- 
damped universality class even when most vibrational 
modes are underdamped. Only a small amount of damp- 
ing is needed to prevent inertia from carrying systems 
over sequential energy barriers, implying that the differ- 
ence between energy barriers is small or the path between 
them complex. Below a critical damping rate a new uni- 
versality class corresponding to the underdamped limit 
is identified. The exponent describing avalanches is close 
to the Gutenberg-Richter law and the finite-size scaling 
function has an unusual form with a plateau before the 
cutoff at large events. Different scaling exponents are ob- 
served at the critical damping rate, indicating that it is 
a multicritical point. 



The scaling exponents in all three regimes (Table QJ 
satisfy the scaling relations /3 + 2a = d and 7 = j3 + 
olt. The hyperscaling relation 7 = d is violated in all 
cases. The number of avalanches at a given energy rises 
less rapidly than system size (7 < d), indicating that 
small events are suppressed by the larger events in bigger 
systems. Exponents obtained from finite-size scaling of 
the distribution of energy and stress drops are consistent. 
However, there is a long power law tail in R(E, L) at 
small E with a size dependent exponent and system-size 
independent cutoff. This tail appears to have dominated 
previous determinations of a llEEi andr [H,^ using 
smaller systems. 

We thank Karin Dahmen for useful discussions. This 
work was supported by the National Science Foundation 
(NSF) under grants DMR-10046442, CMMI-0923018, 
and OCI-108849. 



[1] 
[2] 
[3] 
[4] 
[5] 
[6] 

[7 

[8] 
[9] 

[10 

[11 

[12 

[13 
[14; 
[15 
[16 
[17 
[18 
[19 
[20 



J. Sethna, K. Dahmen, and C. Myers, Nature 410, 242 
(2001), ISSN 0028-0836, 10.1038/35065675. 

B. Gutenberg and C. F. Richter, Bulletin of the Seismo- 
logical Society of America 34, 185 (1944). 

P. A. Houle and J. P. Sethna, Phys. Rev. E 54, 278 
(1996). 

E. Vives, I. Rafols, L. Manosa, J. Ortm, and A. Planes, 
Phys. Rev. B 52, 12644 (1995). 

S. Brinckmann, J.-Y. Kim, and J. R. Greer, Phys. Rev. 
Lett. 100, 155502 (2008). 

M. C. Miguel, A. Vespignani, S. Zapperi, J. Weiss, and 
J.-R. Grasso, Nature 410, 667 (2001), ISSN 0028-0836, 
10.1038/35070524. 

N. Martys, M. O. Robbins, and M. Cieplak, Phys. Rev. 
B 44, 12294 (1991). 

H. Ji and M. O. Robbins, Phys. Rev. B 46, 14519 (1992). 

C. P. C. Prado and Z. Olami, Phys. Rev. A 45, 665 
(1992). 

G. A. Held, D. H. Solina, H. Solina, D. T. Keane, W. J. 
Haag, P. M. Horn, and G. Grinstein, Phys. Rev. Lett. 
65, 1120 (1990). 

H. M. Jaeger, C.-h. Liu, and S. R. Nagel, Phys. Rev. 
Lett. 62, 40 (1989). 

J. M. Carlson and J. S. Langer, Phys. Rev. Lett. 62, 2632 
(1989). 

C. Maloney and A. Lemaitre, Phys. Rev. Lett. 93, 016001 
(2004). 

N. P. Bailey, J. Schi0tz, A. Lemaitre, and K. W. Jacob- 
sen, Phys. Rev. Lett. 98, 095501 (2007). 
H. G. E. Hentschel, S. Karmakar, E. Lerner, and I. Pro- 
caccia, Phys. Rev. Lett. 104, 025501 (2010). 
A. Lemaitre and C. Caroli, Phys. Rev. E 76, 036104 
(2007). 

H. G. E. Hentschel, S. Karmakar, E. Lerner, and I. Pro- 

caccia, Phys. Rev. Lett. 104, 025501 (2010). 

C. E. Maloney and M. O. Robbins, Journal of Physics: 

Condensed Matter 20, 244128 (2008). 

E. Lerner and I. Procaccia, Phys. Rev. E 79, 066109 

(2009). 

V. Privman, Finite size scaling and numerical Simula- 



■5 



tion of statistical systems (World Scientific, 1990), ISBN 
9789810201081. 

[21] L. Landau, E. Lifshitz, A. Kosevich, and L. Pitaevskii, 

Theory of Elasticity, Theoretical Physics (Butterworth- 

Heinemann, 1986), ISBN 9780750626330. 
[22] P. J. Hoogerbrugge and J. M. V. A. Koelman, EPL (Eu- 

rophysics Letters) 19, 155 (1992). 
[23] S. Tewari, D. Schiemann, D. J. Durian, C. M. Knobler, 

S. A. Langer, and A. J. Liu, Phys. Rev. E 60, 4385 

(1999). 



[24] T. Hatano, Phys. Rev. E 79, 050301 (2009). 

[25] K. Dahmen, Y. Ben-Zion, and J. T. Uhl, Nature Physics 
7, 554 (2011), ISSN 1745-2473. 

[26] C. Scholz, The mechanics of earthquakes and fault- 
ing (Cambridge University Press, 2002), ISBN 
9780521655408. 

[27] C. E. Maloney and A. Lemaitre, Phys. Rev. E 74, 016118 
(2006). 



